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A Fast Impurity Solver Based on Gutzwiller variational approach 
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Abstract 

A fast impurity solver for the dynamical mean field theory(DMFT) named Two Mode Approxi- 
mation (TMA) is proposed based on the Gutzwiller variational approach, which captures the main 
features of both the coherent and incoherent motion of the electrons. The new solver works with 
real frequency at zero temperature and it provides directly the spectral function of the electrons. 
It can be easily generalized to multi-orbital impurity problems with general on-site interactions, 
which makes it very useful in LDA+DMFT. Benchmarks on one and two band Hubbard models 
are presented, and the results agree well with those of Exact Diagonalization (ED). 

PACS numbers: 



1 



I. INTRODUCTION 



The accurate calculation of the electronic structure of materials starting from first princi- 
ples is a challenging problem in condensed matter science. The local density approximation 
(LDA) based on density functional theory (DFT) is a widely used ab initio method [1] , which 
has been successfully applied to study the properties of simple metals and semiconductors 
as well as the band insulators. However, it can not be applied to those materials containing 
partially filled narrow bands from d or f shells, because of the so called strong correlation 
effect. 

In LDA the wave like nature rather than the atomic feature of the electronic state is 
emphasized, so it is more suitable to describe those wide energy bands contributed by the 
electrons from outer shells. While for the electrons from those unclosed inner shells hke 3c? 
or 5/ shells, some atomic features such as the multiplet structure remain, which are poorly 
described by LDA. Therefore for those strongly correlated materials, we have to implement 
LDA with some many-body techniques which can deal with the strong correlation effect and 
capture most of the atomic features. 

One notable example of the first- principle schemes is the LDA-I-U method [2] , which can 
successfully describe many interesting effects such as spin, orbital and charge ordering in 
transition metal compounds [3]. Although LDA+U can capture the static orbital and spin 
dependent physics quite well, it still can not consider the dynamical correlation effect, which 
causes lots of interesting phenomena hke Mott transition [4] [5] [6] . 

Another attempt is to use Gutzwiller variational approach [7] [8] to take into account the 
correlation effect (LDA-I-G), which is superior to LDA-I-U and has been successfully applied 
to many systems[9] [10] [11]. LDA-I-G treatment has its advantage in describing ground 
state and low energy excited states, but it can not properly describe the finite temperature 
and dynamical properties due to the lack of high energy excited states. In order to capture 
the overall features of a correlated materials, more sophisticated approaches are needed. 

During the past twenty years, the dynamical mean field theory (DMFT) [12] has been 
quickly developed to be a powerful method to solve the strongly correlated models on the 
lattice. DMFT maps the lattice models to the corresponding quantum impurity models sub- 
ject to self-consistency conditions. Unlike the normal static mean field approaches, DMFT 
keeps the full local dynamics induced by the local interaction. DMFT has been successfully 
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applied to various of correlation problems, such as the Mott transition in Hubbard model 
[13] [14], the pseudo gap behavior in high Tc cuperates [15] and the heavy fermion system 
[16] [17]. Since DMFT can capture quite accurately the correlation feature induced by the 
on-site Coulomb interaction and LDA can take care of the periodic potential as well as the 
long range part of the Coulomb interaction, the combination of the two methods should 
be a very useful scheme for the first priciple calculation of correlation materials. In the 
past twenty years, LDA+DMFT has been developed very quickly and successfully apphed 
to many systems[18], see [19] [20] [21] and [22] for reviews of the recent developments and 
applications. 

In LDA+DMFT, one encounters the problem of how to efficiently solve quantum impurity 
problems with self- consistently determined bath degrees of freedom. A fast impurity solver 
can be regarded as the engine of DMFT, which determines the efficiency and accuracy 
of DMFT. Many impurity solvers have been developed in the past twenty years, which 
can be divided into analytical methods and numerical methods. The analytical mthods 
include equation of motion (EOM) method [23], Hubbard-I approximation [24] [25], iterative 
perturbation theory (IPT) [26] [27], the Non-crossing approximation(NCA) [28] and the 
fluctuation exchange approximation(FLEX) [29]. And the numerical methods include exact 
diagonalization (ED) [30] , Hirsch-Fye Quantum Monte Carlo methods [31] [32] and the 
numerical renormalization group (NRG) [33]. Most recently a powerful continuous-time 
quantum Monte Carlo (CTQMC) solver [34] [35] has also been developed and apphed to 
several reahstic matterials[36] [37]. 

All these impurity solvers have their own advantages and the hmitations as well. Since 
most of the novel quantum phenomena in condensed matter physics happen in very low 
temperature, it is always very important for us to study the low temperature properties 
of the correlated materials using LDA-I-DMFT. Up to now, the impurity solvers which can 
work at extremely low temperature are ED, IPT and NRG. Among them, IPT can only 
apply to the single band system, ED and NRG are numerically quite heavy for a general 
multi-band system. Therefore it is very useful to develop an impurity solver working at zero 
temperature, which satisfles the following criteria, i) It can capture both the low energy 
quasi-particle physics and the high energy Hubbard bands, ii) It works with real frequency 
and gives the real time dynamical properties directly, iii) It is easy to be generalized to 
realistic multi-band systems. 
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Here we propose a fast impurity solver based on Gutzwiller variational approach [9] 
which has the above three advantages. Gutzwiller variational wave function associated 
with Gutzwiller approximation was first proposed to deal with lattice problems such as the 
Hubbard model and the periodical Anderson model [38] [39]. In the present paper, we apply 
a generalized Gutzwiller method called Two Mode Approximation (TMA) to calculate the 
Green's function for a quantum impurity model generated by DMFT. TMA is first proposed 
in reference[40] to calculate the spectral function for the lattice mode. Here we generalize it 
to the quantum impurity problem and make it a useful impurity solver for DMFT. 

In TMA three different types of variational wave functions are constructed for the ground 
states, low energy quasi-particle states and high energy excited state respectively. All the 
variational parameters appearing in different wave functions are determined by minimizing 
the ground state energy, based on which we can obtain the electronic spectral functions over 
the full frequency range. The computational time is mainly determined by the minimization 
of the ground state energy and is similar with the previous study on lattice problem [41], 
which can be easily done even on a single PC. This makes the present approach a fast general 
solver for LDA+DMFT studies. 

The paper is organized as follows. In section II we give the derivation of the method 
and prove that the sum rule for the electronic spectral function is satisfied. In Section HI 
we benchmark our new impurity solver on the two-band Hubbard model with DMFT+ED. 
Finally a summary and the conclusions are made in section IV. 

II. DERIVATION OF THE METHOD 
A. Gutzwiller ground state 

Let us first consider the following multi-orbital impurity Hamiltonian 
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Himp — Hband + Hiocal + Hy 
ka 

Hiocal — ^(Ta'^faTlfa' + SgTlfa- 

<T,a' a 

Hv = J2^ka{dtJa + h.C.) 

ka 

where k denotes the energy levels in the bath and a is the joint index for orbital and spin. In 
Gutzwiller variational approach, the ground state of the above Hamiltonian can be written 
as 

|^) = P|0) (1) 

Where P is the Gutzwiller projector and |0) is a single Slatter Determinant like wave func- 
tion. Both of P and |0) will be determined by minimizing the ground state energy. Following 
reference [9] , the Gutzwiller projector can be written in terms of the projection operators of 
the atomic eigen states as 

P = E^Ar (2) 

In equation (2), the operator jfir = |r)(r| is the projector to the eigen states |r) of 
the atomic Hamiltonian Hiocai, ^-nd are the variational parameters introduced in the 
Gutzwiller theory. Note that if Hiocai only contains density-density interactions, the atomic 
eigen states are known as the Fock states as the following [9], 

r e {0;(l),...,(27V);(l,2),(2,3),...(2iV-l,27V) 

;...(l,..,2Ar)} (3) 

,where N is the number of orbitals. mp is defined as 



= (0|mr|0) (4) 
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Using the operator equalities 

mr = YlnfaY[{l-nfa) (5) 

o-er aeT 

^fo = Xl^r (6) 

with the definition = (0|n/o-|0) and nf„ = (^|n/(^|^), one can prove that — 

n ^/cr n (1 ~ ''^/(t)' ^% = Zlrgo-^^^r = Sr9o-"^r- We would emphasize that 

o-er o-er 

= Ufcr for Gutzwillcr type wave functions with pure density- density interaction, which 
greatly simplify the computation [9, 10]. Therefore the Gutzwiller ground state energy of this 
impurity model reads 

^ ^ (o|m^,p|o) 

' (0|P2|o) ^ ' 

the denominator can be expressed as 

(0|P2|0) = J]mr = l 
r 

while the numerator can be calculated by decomposing the projectors as in equation (5) and 
applying the Wick's theorem [42]. Finally we obtain the ground state energy as 



Eg = J]efe,(0|c+Cfe,|0) + ^£;rmr 

ka r 
+ l];^a^fe.(0|c+/, + /i.C.|0) 



ka 

with 



r9a,r'=rv •\/^/<x(l " ^/a) 
The ground state wave function |^) can be obtained by minimizing the above energy func- 
tional respect to the mr and non-interacting wave function |0)[9, 10] along with the following 
constraints. 

^mr = 1 (8) 
r 

^/c^ = (9) 
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B. zero-temperature Green's function 



For the impurity Hamiltonian Eq.(l), the retarded Green's function for the electrons on 
the impurity site reads 



u + ii] - En + E, 
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^ uj + ir] + Em- Eg ^ ' 

where |^) is the ground state of Himp with the eigen energy Eg, \n) (|m))are the eigen- 
states of Himp with one more (less) electron than the ground state. En and Em arc the 
corresponding eigenvalues. The above expression is exact if the summation of n and m 
includes all the eigenstates. In the present paper, we apply the two mode approximation 
(TMA) to solve the quantum impurity problem, in which we limit the above summation in 
a truncated Hilbert space formed by finite number of excited states over the Gutzwiller vari- 
ational ground state [40, 43] . In order to capture the basic feature of the electronic spectral 
function efficiently, we have to include two types of excited states in TMA, namely the quasi- 
particle excitations which give the right Fermi liquid behavior in low energy, and the high 
energy excited states which are responsible for the Hubbard bands or the atomic multiplet 
features. The former are called quasi-particle states and the latter are called bare-particle 
states in the present paper [40]. The ansatz for the excited states are the following, 

l+M = cLp|0) 
\UHB) = flP\0) 
\QE) = P/t|0) 
\-ka) = Cfe^P|0) 
\LHB) = IP\0) 
\QH) = PU\0) 

where \QE) {\QH)) are the quasi-particle (quasi- hole) states, \UHB) {\LHB) ) are the 
bare-particle (bare-hole) states, and | H — ka) represent the excitations in the bath. 

The excited states listed above are neither orthogonal nor normalized, thus we have to 
calculate the overlaps Oap = {ol\(3) and the matrix elements of the Hamiltonian Ti^/? = 
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{a\H\l3) in this truncated Hilbert space. This procedure could be easily done by applying 
Wick's theorem. We list all the necessary matrix elements and overlaps in the Appendix. 

In order to evaluate the Green's function using expression (10), we have to first obtain 
the cigcn states and eigen values by solving the following generalized eigen equation in the 
truncated Hilbert space. 

n\i) = EiO\i) 

Therefore \l) form a complete basis for the truncated Hilbert space and the completeness 
condition |/)(/| = 1 is satisfied within the truncated Hilbert space. Since both the states 
/^P|0) and /o-P|0) are fully included in the contained Hilbert space, it is easy to prove that 

(--)im[Gr^(^+zry)] = mUl+fium 

TT 

= 1 

, which is the sum rule of the impurity Green's function. 

III. BENCHMARK 

A. Impurity Spectral function 

First of all we check the spectral function obtained by TMA for a single orbital impurity 
model with particle-hole symmetry. The density of states for the heat bath is chosen to be 
the semicircle with the half-width D = 1. The spectral functions for the electron on the 
impurity site with different Hubbard interaction U are shown in Fig.(l). 

Prom Fig.(l) we find that the spectral function contains three parts, the quasi-particle 
peak and two Hubbard bands. With the increment of U, the spectral weight transfers 
from the low energy quasi-particle part to the Hubbard bands. And in large U limit, the 
distance between two Hubbard bands approaches U. All these features are consistent with 
the previous studies on the symmetric Anderson model [44]. In Fig. (2), we compare one 
spectral function for an Anderson impurity model obtained by TMA with that by the normal 
Gutzwiller Approximation (GA)[7, 9] for the lattice model, which only contains the quasi- 
particle part as 
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Compared with normal Gutzwiller approximation (GA lattice), it is very clear that TMA 
can reproduce very nicely the low energy quasi-particle part with slightly smaller spectral 
weight. Therefore the current solver can be viewed as the normal Gutzwiller approximation 
implemented with the Hubbard bands in the high energy part of the electronic spectral 
functions describing the atomic features. 

B. Used as the impurity solver in DMFT 

The present impurity solver can be used in the dynamical mean field theory to study the 
lattice models. In this paper we have studied both the single-band and two-band Hubbard 
model at paramagnetic phase with arbitrary fillings. 

1. Single-hand Hubbard model 

Wc start with the single band Hubbard model on the Bethe lattice with half band width 
D — 1. First we check the half filling case. We show the spectral function with the increment 
of U in Fig. (3), from which we see that the height of quasi-particle peak changes httle before 
Mott transition, but the integral of the quasi-particle spectrum reduces as U increases. This 
feature is consistent with the previous results obtained by DMFT-|-IPT[12]. 

We show the results for the systems away from half filling in Fig. (4). 

With the increment of filling factor from Ntot = 0.2 to half filling Ntot = 1-0, the spectral 
weight continuously transfers from the low energy quasi-particle part to the high energy Hub- 
bard bands, which is consistent with the common understanding that the strong correlation 
effect is less pronounced when the system is doped away from half filling. 

In Fig. (5), we quantitatively compare the density of states (DOS) obtained by 
DMFT-FTMA with that by DMFT-FED. We find quite good agreement between them for 
both the half filling and non-half filling cases. While we also find two disagreements. Com- 
pared with the DMFT-I-ED results, the total spectral weight of the quasi-particle part is 
over-estimated while the width of the Hubbard bands is under-estimated by DMFT-I-TMA. 
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We have also calculated the quasi-particle weight which is a characteristic quantity 
describing the strength of the correlation effect and is defined as: 

dRe[E^{uj + iri)] 

^0- = (1 ) |a;=0 (l^j 

In Fig. (6) we show quasi-particle weight obtained by DMFT+TMA as the function of U 
for different fiUing factors. In the half filling case, the value of z decreases as the increment 
of U until the critical Uc for the Mott transition. As shown in Fig. (6), Uc obtained by 
DMFT+TMA is around 3.6,which is bigger than Uc2 = 2.9 obtained by DMFT+ED. 

In Fig. (7), we compare the 2;- factors obtained by DMFT+TMA, Gutzwiller approxima- 
tion on the lattice model (lattice GA) and DMFT+ED. As discussed in reference[10] and 
[40] , we can only obtain the ground state energy quite accurately by lattice GA, but not for 
the 2;-factor. The reason is quite obvious that in the lattice GA only the low energy quasi- 
particle states in equation(lO) can be considered, which limits the accuracy of 2;-factor. 
While in TMA, we first apply the DMFT scheme to treat the inter-site correlation on a 
mean field level, which is in principle similar with GA. Then in solving the effective impu- 
rity model, we enlarge the variational space by including more excited states, which gives 
us more accurate description of the low energy excited states and reduces the disagreement 
in 2;-factor with DMFT+ED results as shown in Fig. (7). 

2. Two-band Hubbard model on the Bethe lattice 

The situation becomes more complicated when we consider two-band models. We start 
with the simplest case that the two bands are degenerate with half bandwidth Di = D2 — 1 
and the local part of the Hamiltonian has SU(4) symmetry, which can be written as 

Hat ^ U'^hb,^hb^l + U'^hi^^h2,a' (13) 

b <7,a' 

We first show the quasi-particle weight obtained by DMFT+TMA versus U at different 
filling factors in Fig(8) and the comparison with DMFT+ED and lattice GA in Fig(9). 

The Mott transition at integer fillings can be observed with Uc slightly larger than the 
DMFT+ED results. As shown in Fig(9), the improvement of the quasi-particle weight 
against the lattice GA is quite dramatic, which indicates that even for the low energy quasi- 
particle part the DMFT+TMA is better than applying the GA directly to the lattice model. 
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The behavior of z as the function of the fiUing factor for fixed IJ = 5.0 is shown in 
Fig(lO), from which we can find that compared with lattice GA the results obtained by 
DMFT+TMA is much closer to DMFT+ED. 

Next we take the Hund's coupling constant J into account. Then the atomic Hamiltonian 
becomes 

Hat = U^Y^ nb,^nb,l + U'^ ni^afl2,a' - fli^afl2,a 
b a,(7' (T 

a 

+c+^c+^Ci4Ci,t) (14) 

We have the relation U — U' = 2J for system with cubic symmetry[45]. In the current 
study, we only keep the longitudinal part of the Hund's rule coupling and neglect the spin flip 
and pair hopping terms which correspond to the last two terms in the above equation. The 
results for the full rotational invariance interaction will be studied in detail and published 
elsewhere. 

The quasi-particle weight obtained by DMFT+TMA as the function of U is shown in 
Fig. (11). We also compare the results with DMFT+ED in Fig. (12), from which we find that 
Uc obtained from TMA is larger than that of DMFT+ED as for the single band model. 

In Fig. (11), we find that the Brinkman-Rice(BR) transition is continuous only at the point 
Jz — and first order like for all non-zero J^, which is similar with the results in reference 
[9] obtained by rotational invariant Gutzwiller approximation. This similarity indicates 
that for degenerate multi-band Hubbard model the basic feature of the BR transition docs 
not strongly relics on the variational invariant treatment of the interaction. Moreover, the 
similar discontinuity and the tendency that the critical Uc decreases as Jz/U increases is 
also obtained in [46] , where the self-energy functional method is used. 

However, for the non-degenerate multi-band models, i.e. the two-band model with differ- 
ent band widths, the correct variational invariant treatment is necessary to obtain some of 
the qualitative features like the orbital selective Mott transition (OSMT)[47]. The detailed 
study for the OSMT using the variational invariant TMA solver will be presented elsewhere. 
Here we only give the results for an extreme case, where the band width difference of the two 
bands is very large. In Fig. (13) and (14), we represent the DOS as well as the quasi-particle 
weight as the functional of U with fixed Jz/U — 0.3 and half band width Di — 1.0, D2 — 6.0. 
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Obviously in such extreme case, the system is in the orbital selective Mott phase which is 
consistent with reference [48]. 



IV. CONCLUSIONS 

In this paper wc present a new impurity solver named Two Mode Approximation (TMA) 
for the multi-orbital quantum impurity model generated by DMFT. By constructing the 
trial wave functions based on the Gutzwiller variational theory not only for the ground state 
but also the low energy and high energy excited states, we can obtain the spectral functions 
of the electrons on the impurity level with the satisfactory of the sum rule. Compared with 
other popular impurity solvers, TMA works with the real frequency and can obtain both 
the low energy quasi-particle and high energy Hubbard band behavior. Moreover TMA can 
be generalized to treat the problem with quite general on-site interaction, which make it a 
good solver to be used in LDA-^DMFT. 

ACKNOWLEDGEMENT: The authors would thank Q. M. Liu, X. Y. Deng, Y. Wan and 
N.H. Tong for their helpful discussions. We acknowledge the supports from NSF of China , 
and that from the 973 program of China (No.2007CB925000). 

V. APPENDIX: OVERLAPS AND HAMILTONIAN ELEMENTS 
A. Overlaps 



Define 




the non- vanishing overlaps are 



{+kia \+k2(T) 



mk^a^t, |0) 



{+ka\UHB) 



z,{0\ckJl |0) 



i+ka \QE) 



^kJl |0) 
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{UHB\UHB) = {l-n%) 
{UHB\QE)^z,{l-n%) 
{QE\QE) = {l-n%) 
{-kia\-k2a) = (Olc+^Cfe^a |0) 
{-ka\LHB) = z,{0\cU^\0) 
{-ka\QH)^{0\ctS\0) 
{LHB\LHB) = n}^ 
{LHB\QH) = z,n% 
{LHB\LHB) = n% 

B. Hamiltonian Elements 



H 

Hband 
Hiocal 

Hv 



Hband + Hiocal + Hy 
ka 

r <T rgo- 



ka 
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1. H_band 



Define 



and 



XfTtr' 



j2 - 



_ v^"'ri/"i-, 

r,3.,r,5.' Jn%{l-n%)n%,{l-n%,) 
r2=riU(T'\o- ^ 



r25.,r.5.' Jn%{l - n%)n%,{l - n%,) 
ri=r2Uo-uo-' * 

ri3a,ri9<^ (1 - i^%)\/n%,{l-nj^,) 
r2=riua' J y J 



era' / J 



we will have 



V3a,V3u' /'^ 



^0 n - ) 



(+A;i(7 



band 



+ k2(T 



efcV'[5a<x'(0|Cfc,^C+^Cfc/^C+^ |0) 

k'a' 

+ (1 - (^aaO 

x(5+t(0|c,,.c+ Jt/. |0) (0|c+,,cfcv'/i/.' |0) 

+i?+,7(0|cfe,.c+ Jt/^ |0) (o|c+,,c,v'/.'/i |0) 
+S+;(0|c,,.c+ Jjt |o) (0|c+,,c,,.,/J,/., |0) 

+B;;;{o\c,,,ctJJl |o) (o|c+,,Cfc,., |o))] 
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+ (1 - 5,,,){0\CkJl |0) {xMOlK'a'^k'a'fi'fa' |0> + ^;av(0|c+^,CfeV'/a'/J. |0))] 



+(1 - S,,,){0\ckJl |0) (S++(0|c+,,Cfev'/J,/.' |0) + S+7(0|c+,,c,v'/.'& |0))] 



+(i-5,,0(i-O 

x(S+;(0|c+,,c,v'/i/.' |0) +S;,7(0|c+,,Cfev'/.'/J. |0)) 



+(l-(5,,0(l-O 

x(a;<,v(0|c+^,cfcv'/i/a' |0) + v<,v(0|c+^,Cfev//^//i |0)) 



{QE\Hband\QE) = Y.'k'ASaa'{0\fadt'aCk'afUO) 



x(B:^,{o\ct,Ak'a'flU' |0) + B:;,{o\ct,Ak'a'U'fi' m 



+(1 - 

ABt:>mtAJlh |0) (0|c+,,c,v'/i/.' |0) 

+i?:,7(0|c+,Cfe, Jt/^ |0) (0|c+,,c,v'/.'/i |0) 

+ByM^lAkJJl |0) (0|c+,,c,v'/i/.' |0) 
+B-;,{Q\clA,JJl |0) (0|c+,,c,v'/.'/J, |0»] 
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kcr Hband\LHB) = ^ek'a'{Sc7a'Za{0\CkaCt'a^k'af<j 



k'a' 



+ (1 - S„^,){0\C+J^ |0) {x^'AO\ct'a'Ck'a'fUa' |0) + V„,„{0\c+^,Ck>., |0))) 



—ka 



k'a' 



+ (1 - S,,>){0\cll |0) (i?J(0|c+,,c,v'/i/.' |0) + i?,,7(0|c+,,c,v'/.'/i |0))) 



{LHB\Hba„d\LHB) = efev'[<^aa'(0|/]a^vCfev/a |0) 

+(1 - 6,^,)n% 

x{Bt:>m,,.c,.,'flh' |0) + i^+7(0|c+,,CfeV'/.'/i |0)) 



k'a' 

+(1 - S^^>)n% 

x{Xa'a{0\c'^,^,Ck'a'fl'fa' |0) + Va'a{0\ct'^'Ck'a' fa' ft' |0)) 



{QH\Hband\QH) = J2'k'a'[Saa'{0\flct,^Ck'afa\0) 
k'a' 

+ (1 - 

x(i?,-+(0|c,V4v'/;/.' |0) + B^;,{0\ct,,,c,,,'fa'fl, |0)) 



5. HJocal 

Here we define a function for set: 

1, (7 e r 

Aa,r = { 
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a' 



Si = ^ Ernir + ^ ^ mr 

r TBa' 



r 

= y^^A^,^^yrn^m^\aS^ 



r 



r 

= y^ A^,rmr\a>S'r 



7) = y^ ■£;r(l - ^a,r) ymrmrua + y^ y^(l - ^a,r) Vmrmrug 
r 

= y^(l - ^a,r)v'mrmrutT'5'r 



r 



r 
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+ E E A'r',r[Aa,r^ - (1 - Aa,r)j^ 
= K.r;^-(l-^..r)r^lSr 



thus 



{+ha Hiocall + k2a) = {0\Ck,aC+ |0) S, + {0\Ck,Ji |0) {0\Uc+ |0) ,S25((7) 



Hlocal\UHB 



{+ka HiocailQE) = —{0\ckan\0)Ss{a) 



n 



{UHB\Hiocai\UHB)^S4{a) 



{UHB\Hiocai\QE) 
{QE\Hi,^i\QE) 



1-71°. 



n 



{-ha Hiocail - ha) = {0\cl„Ck,a |0) Si - (0|c+ J, |0) {0\flck,. |0) 525(a) 



{-ka 



Hlocal\LHB ) — 



■■mja\0)SM 



(-ka 



HlocallQH 



\ctJa\0)S,{a) 



{LHB\Hiocai\LHB) = Sj{a) 



{LHB\Hi,,,i\QH) 



n 



-S.ia) 
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{QH\Hiocai\QH) = -^Se{<j) 



3. H_V 



{+kia\Hv\ + k2(T) = ^Vk'a'[5aa'Za{{^\CkiaCl,JaCt^a |0) + (OlCfeia/lcfeVCj^^a |0)) 

fc'cr' 

+ (1 - 5.,0((0|4v/.' |0) + (0|/J,c,v' |0)) 
x(x,,,(0|c,,,c+ Jt/^ |0) + v,,,{Q\ck,.ctJJl |0))] 



{-^ka Hv\UHBj = J2Vk',,[6,.,{0\ckactJJl\0) 

k'a' 

+ (1 - 6^^,){0\ckJi |0) (w(0|c+,,A. |0) + w,,>{0\fl,Ck'a' |0))] 



k'a' 

+ + (1 - 5..0w(o|cfe./i |o) m4,,u |o) + (o|/J,CfeV' |o))] 



{UHB\Hv\UHB) = J2 yk'a'il - S.a')Vaa'{l ' O((0|c+,,A, |0) + (Ol/^Cfev |0)) 



k'a' 



{UHB\Hv\QE) = J2 Vk'a'il - Saa')(l - nJJK,,(0|c+^,/,, |0) + yaa'{0\flck'a' m 



k'a' 



{QE\Hy\QE) = ^^'-'(1 - '^-Owll - n%)m4^,U, |0) + (0|/],Cfev' |0)) 
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(-/CiCr Hv\-k2aj = ^Vk'a'[5aa'Za{{0\ct^M,JaCk2a\0) + {0\ct,jlck'aC^^^ 

k'a' 

+(1 - 5..0((o|a^v'/<x' |o) + (o|/;,Cfc,., |o)) 

{x,,,{0\cl,Ck,M, |0) + w(0|c+,Cfe, JJt |o))] 



[-ka 



Hy\LHB) = 5^VfcvKa'(0|c+Jta^v/cr|0) 
k'a' 

+(1 - 5..0(0|5t/. |0) K..(0|c+,,/., |0) + y,..(0|/;cfev' |0))] 



fc'tr' 

+ (1 - <^..')w(0|c+ |0) ((0|4,,/,, |0) + (0|/J,Cfev' |0))] 



{LHB\Hv\LHB) = - 5aa')w<((0|c+,,/,, |0) + (0|/icfcv' |0» 



k'a' 



{LHB\Hv\QH) = - <^.a')^/a(y.a'(0|4v'^' |0) + «'..'(0|/icfev' |0)) 



fe'tr' 



(g// \Hv\QH) = J] - <^aa')wn;,((0|c+,,/,, |0) + (0|/J:,Cfcv' |0)) 



k'a' 
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FIG. 1: The spectral function of electrons on the impurity site for an single orbital impurity model 
with different U and semi-circular density of states in the bath. 



FIG. 2: The spectral function of electrons on the impurity site for an single orbital impurity model 
obtained by TMA and GA lattice with U = 1. 
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FIG. 3: The density of states (DOS) obtained by DMFT+TMA for single-band Hubbard model 
on Bethe lattice at half filling. 
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FIG. 4: The density of states (DOS) obtained by DMFT+TMA of single band Hubbard model 
under U = 2 with different fillings. 
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FIG. 6: Quasi-particle weight z of single band Hubbard model obtained by DMFT+TMA versus 
U at different fillings. 
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FIG. 7: Comparison of Quasi-particle weight z for the single band Hubbard model obtained by 
DMFT+TMA, GA lattice and DMFT+ED. 
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FIG. 8: Quasi-particle weight z as the function of U for the two-band Hubbard model with SU(N) 
symmetry obtained by DMFT-I-TMA. 
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FIG. 9: Comparison of quasi-particle weight z for the two band Hubbard model with SU(N) 
symmetry at different filhngs obtained by DMFT+TMA, GA lattice and DMFT+ED. 



— GAJattice 
-•-TMA 
-i- DMFTED 
U=5, twoband SU(N) 




FIG. 10: Comparison of quasi-particle weight z as the function of total number of particles for the 
two-band Hubbard model with SU(N) symmetry at [/ = 5 obtained by DMFT-FTMA, GA lattice 
and DMFT-hED. 
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FIG. 11: Quasi-particle weight z as the function of U for the two-degenerate-band Hubbard model 
with longitudinal Hund's coupling obtained by DMFT-I-TMA at different Jz/U . 
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FIG. 12: Comparison of quasi-particle weight z for the two-degenerate-band Hubbard model with 
longitudinal Hund's couphng obtained by DMFT-FTMA and DMFT-hED at different Jz/U. 
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FIG. 13: The spectral functions obtained by DMFT+TMA for two-nondegenerate-band Hubbard 
model with band width ratio 1 : 6 under different U with = 0.3U. 
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FIG. 14: Quasi-particle weight z of different bands as the function of U for two-nondegenerate-band 
Hubbard model with band width ratio 1 : 6 under different U with = 0.3U. 
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A Fast Impurity Solver Based on Gutzwiller variational approach 

Jia Ning Zhuang, Lei Wang, Zhong Fang, Xi Dai 
(Dated: February 8, 2009) 

Abstract 

A fast impurity solver for the dynamical mean field theory(DMFT) named Two Mode Approxi- 
mation (TMA) is proposed based on the Gutzwiller variational approach, which captures the main 
features of both the coherent and incoherent motion of the electrons. The new solver works with 
real frequency at zero temperature and it provides directly the spectral function of the electrons. 
It can be easily generalized to multi-orbital impurity problems with general on-site interactions, 
which makes it very useful in LDA+DMFT. Benchmarks on one and two band Hubbard models 
are presented, and the results agree well with those of Exact Diagonalization (ED). 

PACS numbers: 
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I. INTRODUCTION 



The accurate calculation of the electronic structure of materials starting from first princi- 
ples is a challenging problem in condensed matter science. The local density approximation 
(LDA) based on density functional theory (DFT) is a widely used ab initio method [1] , which 
has been successfully applied to study the properties of simple metals and semiconductors 
as well as the band insulators. However, it can not be applied to those materials containing 
partially filled narrow bands from d or f shells, because of the so called strong correlation 
effect. 

In LDA the wave like nature rather than the atomic feature of the electronic state is 
emphasized, so it is more suitable to describe those wide energy bands contributed by the 
electrons from outer shells. While for the electrons from those unclosed inner shells hke 3c? 
or 5/ shells, some atomic features such as the multiplet structure remain, which are poorly 
described by LDA. Therefore for those strongly correlated materials, we have to implement 
LDA with some many-body techniques which can deal with the strong correlation effect and 
capture most of the atomic features. 

One notable example of the first- principle schemes is the LDA-I-U method [2] , which can 
successfully describe many interesting effects such as spin, orbital and charge ordering in 
transition metal compounds [3]. Although LDA+U can capture the static orbital and spin 
dependent physics quite well, it still can not consider the dynamical correlation effect, which 
causes lots of interesting phenomena hke Mott transition [4] [5] [6] . 

Another attempt is to use Gutzwiller variational approach [7] [8] to take into account the 
correlation effect (LDA-I-G), which is superior to LDA-I-U and has been successfully applied 
to many systems[9] [10] [11]. LDA-I-G treatment has its advantage in describing ground 
state and low energy excited states, but it can not properly describe the finite temperature 
and dynamical properties due to the lack of high energy excited states. In order to capture 
the overall features of a correlated materials, more sophisticated approaches are needed. 

During the past twenty years, the dynamical mean field theory (DMFT) [12] has been 
quickly developed to be a powerful method to solve the strongly correlated models on the 
lattice. DMFT maps the lattice models to the corresponding quantum impurity models sub- 
ject to self-consistency conditions. Unlike the normal static mean field approaches, DMFT 
keeps the full local dynamics induced by the local interaction. DMFT has been successfully 
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applied to various of correlation problems, such as the Mott transition in Hubbard model 
[13] [14], the pseudo gap behavior in high Tc cuperates [15] and the heavy fermion system 
[16] [17]. Since DMFT can capture quite accurately the correlation feature induced by the 
on-site Coulomb interaction and LDA can take care of the periodic potential as well as the 
long range part of the Coulomb interaction, the combination of the two methods should 
be a very useful scheme for the first priciple calculation of correlation materials. In the 
past twenty years, LDA+DMFT has been developed very quickly and successfully apphed 
to many systems[18], see [19] [20] [21] and [22] for reviews of the recent developments and 
applications. 

In LDA+DMFT, one encounters the problem of how to efficiently solve quantum impurity 
problems with self- consistently determined bath degrees of freedom. A fast impurity solver 
can be regarded as the engine of DMFT, which determines the efficiency and accuracy 
of DMFT. Many impurity solvers have been developed in the past twenty years, which 
can be divided into analytical methods and numerical methods. The analytical mthods 
include equation of motion (EOM) method [23], Hubbard-I approximation [24] [25], iterative 
perturbation theory (IPT) [26] [27], the Non-crossing approximation(NCA) [28] and the 
fluctuation exchange approximation(FLEX) [29]. And the numerical methods include exact 
diagonalization (ED) [30] , Hirsch-Fye Quantum Monte Carlo methods [31] [32] and the 
numerical renormalization group (NRG) [33]. Most recently a powerful continuous-time 
quantum Monte Carlo (CTQMC) solver [34] [35] has also been developed and apphed to 
several reahstic matterials[36] [37]. 

All these impurity solvers have their own advantages and the hmitations as well. Since 
most of the novel quantum phenomena in condensed matter physics happen in very low 
temperature, it is always very important for us to study the low temperature properties 
of the correlated materials using LDA-I-DMFT. Up to now, the impurity solvers which can 
work at extremely low temperature are ED, IPT and NRG. Among them, IPT can only 
apply to the single band system, ED and NRG are numerically quite heavy for a general 
multi-band system. Therefore it is very useful to develop an impurity solver working at zero 
temperature, which satisfles the following criteria, i) It can capture both the low energy 
quasi-particle physics and the high energy Hubbard bands, ii) It works with real frequency 
and gives the real time dynamical properties directly, iii) It is easy to be generalized to 
realistic multi-band systems. 
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Here we propose a fast impurity solver based on Gutzwiller variational approach [9] 
which has the above three advantages. Gutzwiller variational wave function associated 
with Gutzwiller approximation was first proposed to deal with lattice problems such as the 
Hubbard model and the periodical Anderson model [38] [39]. In the present paper, we apply 
a generalized Gutzwiller method called Two Mode Approximation (TMA) to calculate the 
Green's function for a quantum impurity model generated by DMFT. TMA is first proposed 
in reference[40] to calculate the spectral function for the lattice mode. Here we generalize it 
to the quantum impurity problem and make it a useful impurity solver for DMFT. 

In TMA three different types of variational wave functions are constructed for the ground 
states, low energy quasi-particle states and high energy excited state respectively. All the 
variational parameters appearing in different wave functions are determined by minimizing 
the ground state energy, based on which we can obtain the electronic spectral functions over 
the full frequency range. The computational time is mainly determined by the minimization 
of the ground state energy and is similar with the previous study on lattice problem [41], 
which can be easily done even on a single PC. This makes the present approach a fast general 
solver for LDA+DMFT studies. 

The paper is organized as follows. In section II we give the derivation of the method 
and prove that the sum rule for the electronic spectral function is satisfied. In Section HI 
we benchmark our new impurity solver on the two-band Hubbard model with DMFT+ED. 
Finally a summary and the conclusions are made in section IV. 

II. DERIVATION OF THE METHOD 
A. Gutzwiller ground state 

Let us first consider the following multi-orbital impurity Hamiltonian 
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Himp — Hband + Hiocal + Hy 
ka 

Hiocal — ^(Ta'^faTlfa' + SgTlfa- 

<T,a' a 

Hv = J2^ka{dtJa + h.C.) 

ka 

where k denotes the energy levels in the bath and a is the joint index for orbital and spin. In 
Gutzwiller variational approach, the ground state of the above Hamiltonian can be written 
as 

|^) = P|0) (1) 

Where P is the Gutzwiller projector and |0) is a single Slatter Determinant like wave func- 
tion. Both of P and |0) will be determined by minimizing the ground state energy. Following 
reference [9] , the Gutzwiller projector can be written in terms of the projection operators of 
the atomic eigen states as 

P = E^Ar (2) 

In equation (2), the operator jfir = |r)(r| is the projector to the eigen states |r) of 
the atomic Hamiltonian Hiocai, ^-nd are the variational parameters introduced in the 
Gutzwiller theory. Note that if Hiocai only contains density-density interactions, the atomic 
eigen states are known as the Fock states as the following [9], 

r e {0;(l),...,(27V);(l,2),(2,3),...(2iV-l,27V) 

;...(l,..,2Ar)} (3) 

,where N is the number of orbitals. mp is defined as 



= (0|mr|0) (4) 
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Using the operator equalities 

mr = YlnfaY[{l-nfa) (5) 

o-er aeT 

^fo = Xl^r (6) 

with the definition = (0|n/o-|0) and nf„ = (^|n/(^|^), one can prove that — 

n ^/cr n (1 ~ ''^/(t)' ^% = Zlrgo-^^^r = Sr9o-"^r- We would emphasize that 

o-er o-er 

= Ufcr for Gutzwillcr type wave functions with pure density- density interaction, which 
greatly simplify the computation [9, 10]. Therefore the Gutzwiller ground state energy of this 
impurity model reads 

^ ^ (o|m^,p|o) 

' (0|P2|o) ^ ' 

the denominator can be expressed as 

(0|P2|0) = J]mr = l 
r 

while the numerator can be calculated by decomposing the projectors as in equation (5) and 
applying the Wick's theorem [42]. Finally we obtain the ground state energy as 



Eg = J]efe,(0|c+Cfe,|0) + ^£;rmr 

ka r 
+ l];^a^fe.(0|c+/, + /i.C.|0) 



ka 

with 



r9a,r'=rv •\/^/<x(l " ^/a) 
The ground state wave function |^) can be obtained by minimizing the above energy func- 
tional respect to the mr and non-interacting wave function |0)[9, 10] along with the following 
constraints. 

^mr = 1 (8) 
r 

^/c^ = (9) 
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B. zero-temperature Green's function 



For the impurity Hamiltonian Eq.(l), the retarded Green's function for the electrons on 
the impurity site reads 



u + ii] - En + E, 
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^ uj + ir] + Em- Eg ^ ' 

where |^) is the ground state of Himp with the eigen energy Eg, \n) (|m))are the eigen- 
states of Himp with one more (less) electron than the ground state. En and Em arc the 
corresponding eigenvalues. The above expression is exact if the summation of n and m 
includes all the eigenstates. In the present paper, we apply the two mode approximation 
(TMA) to solve the quantum impurity problem, in which we limit the above summation in 
a truncated Hilbert space formed by finite number of excited states over the Gutzwiller vari- 
ational ground state [40, 43] . In order to capture the basic feature of the electronic spectral 
function efficiently, we have to include two types of excited states in TMA, namely the quasi- 
particle excitations which give the right Fermi liquid behavior in low energy, and the high 
energy excited states which are responsible for the Hubbard bands or the atomic multiplet 
features. The former are called quasi-particle states and the latter are called bare-particle 
states in the present paper [40]. The ansatz for the excited states are the following, 

l+M = cLp|0) 
\UHB) = flP\0) 
\QE) = P/t|0) 
\-ka) = Cfe^P|0) 
\LHB) = IP\0) 
\QH) = PU\0) 

where \QE) {\QH)) are the quasi-particle (quasi- hole) states, \UHB) {\LHB) ) are the 
bare-particle (bare-hole) states, and | H — ka) represent the excitations in the bath. 

The excited states listed above are neither orthogonal nor normalized, thus we have to 
calculate the overlaps Oap = {ol\(3) and the matrix elements of the Hamiltonian Ti^/? = 
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{a\H\l3) in this truncated Hilbert space. This procedure could be easily done by applying 
Wick's theorem. We list all the necessary matrix elements and overlaps in the Appendix. 

In order to evaluate the Green's function using expression (10), we have to first obtain 
the cigcn states and eigen values by solving the following generalized eigen equation in the 
truncated Hilbert space. 

n\i) = EiO\i) 

Therefore \l) form a complete basis for the truncated Hilbert space and the completeness 
condition |/)(/| = 1 is satisfied within the truncated Hilbert space. Since both the states 
/^P|0) and /o-P|0) are fully included in the contained Hilbert space, it is easy to prove that 

(--)im[Gr^(^+zry)] = mUl+fium 

TT 

= 1 

, which is the sum rule of the impurity Green's function. 

III. BENCHMARK 

A. Impurity Spectral function 

First of all we check the spectral function obtained by TMA for a single orbital impurity 
model with particle-hole symmetry. The density of states for the heat bath is chosen to be 
the semicircle with the half-width D = 1. The spectral functions for the electron on the 
impurity site with different Hubbard interaction U are shown in Fig.(l). 

Prom Fig.(l) we find that the spectral function contains three parts, the quasi-particle 
peak and two Hubbard bands. With the increment of U, the spectral weight transfers 
from the low energy quasi-particle part to the Hubbard bands. And in large U limit, the 
distance between two Hubbard bands approaches U. All these features are consistent with 
the previous studies on the symmetric Anderson model [44]. In Fig. (2), we compare one 
spectral function for an Anderson impurity model obtained by TMA with that by the normal 
Gutzwiller Approximation (GA)[7, 9] for the lattice model, which only contains the quasi- 
particle part as 
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Compared with normal Gutzwiller approximation (GA lattice), it is very clear that TMA 
can reproduce very nicely the low energy quasi-particle part with slightly smaller spectral 
weight. Therefore the current solver can be viewed as the normal Gutzwiller approximation 
implemented with the Hubbard bands in the high energy part of the electronic spectral 
functions describing the atomic features. 

B. Used as the impurity solver in DMFT 

The present impurity solver can be used in the dynamical mean field theory to study the 
lattice models. In this paper we have studied both the single-band and two-band Hubbard 
model at paramagnetic phase with arbitrary fillings. 

1. Single-hand Hubbard model 

Wc start with the single band Hubbard model on the Bethe lattice with half band width 
D — 1. First we check the half filling case. We show the spectral function with the increment 
of U in Fig. (3), from which we see that the height of quasi-particle peak changes httle before 
Mott transition, but the integral of the quasi-particle spectrum reduces as U increases. This 
feature is consistent with the previous results obtained by DMFT-|-IPT[12]. 

We show the results for the systems away from half filling in Fig. (4). 

With the increment of filling factor from Ntot = 0.2 to half filling Ntot = 1-0, the spectral 
weight continuously transfers from the low energy quasi-particle part to the high energy Hub- 
bard bands, which is consistent with the common understanding that the strong correlation 
effect is less pronounced when the system is doped away from half filling. 

In Fig. (5), we quantitatively compare the density of states (DOS) obtained by 
DMFT-FTMA with that by DMFT-FED. We find quite good agreement between them for 
both the half filling and non-half filling cases. While we also find two disagreements. Com- 
pared with the DMFT-I-ED results, the total spectral weight of the quasi-particle part is 
over-estimated while the width of the Hubbard bands is under-estimated by DMFT-I-TMA. 
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We have also calculated the quasi-particle weight which is a characteristic quantity 
describing the strength of the correlation effect and is defined as: 

dRe[E^{uj + iri)] 

^0- = (1 ) |a;=0 (l^j 

In Fig. (6) we show quasi-particle weight obtained by DMFT+TMA as the function of U 
for different fiUing factors. In the half filling case, the value of z decreases as the increment 
of U until the critical Uc for the Mott transition. As shown in Fig. (6), Uc obtained by 
DMFT+TMA is around 3.6,which is bigger than Uc2 = 2.9 obtained by DMFT+ED. 

In Fig. (7), we compare the 2;- factors obtained by DMFT+TMA, Gutzwiller approxima- 
tion on the lattice model (lattice GA) and DMFT+ED. As discussed in reference[10] and 
[40] , we can only obtain the ground state energy quite accurately by lattice GA, but not for 
the 2;-factor. The reason is quite obvious that in the lattice GA only the low energy quasi- 
particle states in equation(lO) can be considered, which limits the accuracy of 2;-factor. 
While in TMA, we first apply the DMFT scheme to treat the inter-site correlation on a 
mean field level, which is in principle similar with GA. Then in solving the effective impu- 
rity model, we enlarge the variational space by including more excited states, which gives 
us more accurate description of the low energy excited states and reduces the disagreement 
in 2;-factor with DMFT+ED results as shown in Fig. (7). 

2. Two-band Hubbard model on the Bethe lattice 

The situation becomes more complicated when we consider two-band models. We start 
with the simplest case that the two bands are degenerate with half bandwidth Di = D2 — 1 
and the local part of the Hamiltonian has SU(4) symmetry, which can be written as 

Hat ^ U'^hb,^hb^l + U'^hi^^h2,a' (13) 

b <7,a' 

We first show the quasi-particle weight obtained by DMFT+TMA versus U at different 
filling factors in Fig(8) and the comparison with DMFT+ED and lattice GA in Fig(9). 

The Mott transition at integer fillings can be observed with Uc slightly larger than the 
DMFT+ED results. As shown in Fig(9), the improvement of the quasi-particle weight 
against the lattice GA is quite dramatic, which indicates that even for the low energy quasi- 
particle part the DMFT+TMA is better than applying the GA directly to the lattice model. 
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The behavior of z as the function of the fiUing factor for fixed IJ = 5.0 is shown in 
Fig(lO), from which we can find that compared with lattice GA the results obtained by 
DMFT+TMA is much closer to DMFT+ED. 

Next we take the Hund's coupling constant J into account. Then the atomic Hamiltonian 
becomes 

Hat = U^Y^ nb,^nb,l + U'^ ni^afl2,a' - fli^afl2,a 
b a,(7' (T 

a 

+c+^c+^Ci4Ci,t) (14) 

We have the relation U — U' = 2J for system with cubic symmetry[45]. In the current 
study, we only keep the longitudinal part of the Hund's rule coupling and neglect the spin flip 
and pair hopping terms which correspond to the last two terms in the above equation. The 
results for the full rotational invariance interaction will be studied in detail and published 
elsewhere. 

The quasi-particle weight obtained by DMFT+TMA as the function of U is shown in 
Fig. (11). We also compare the results with DMFT+ED in Fig. (12), from which we find that 
Uc obtained from TMA is larger than that of DMFT+ED as for the single band model. 

In Fig. (11), we find that the Brinkman-Rice(BR) transition is continuous only at the point 
Jz — and first order like for all non-zero J^, which is similar with the results in reference 
[9] obtained by rotational invariant Gutzwiller approximation. This similarity indicates 
that for degenerate multi-band Hubbard model the basic feature of the BR transition docs 
not strongly relics on the variational invariant treatment of the interaction. Moreover, the 
similar discontinuity and the tendency that the critical Uc decreases as Jz/U increases is 
also obtained in [46] , where the self-energy functional method is used. 

However, for the non-degenerate multi-band models, i.e. the two-band model with differ- 
ent band widths, the correct variational invariant treatment is necessary to obtain some of 
the qualitative features like the orbital selective Mott transition (OSMT)[47]. The detailed 
study for the OSMT using the variational invariant TMA solver will be presented elsewhere. 
Here we only give the results for an extreme case, where the band width difference of the two 
bands is very large. In Fig. (13) and (14), we represent the DOS as well as the quasi-particle 
weight as the functional of U with fixed Jz/U — 0.3 and half band width Di — 1.0, D2 — 6.0. 
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Obviously in such extreme case, the system is in the orbital selective Mott phase which is 
consistent with reference [48]. 



IV. CONCLUSIONS 

In this paper wc present a new impurity solver named Two Mode Approximation (TMA) 
for the multi-orbital quantum impurity model generated by DMFT. By constructing the 
trial wave functions based on the Gutzwiller variational theory not only for the ground state 
but also the low energy and high energy excited states, we can obtain the spectral functions 
of the electrons on the impurity level with the satisfactory of the sum rule. Compared with 
other popular impurity solvers, TMA works with the real frequency and can obtain both 
the low energy quasi-particle and high energy Hubbard band behavior. Moreover TMA can 
be generalized to treat the problem with quite general on-site interaction, which make it a 
good solver to be used in LDA-^DMFT. 
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V. APPENDIX: OVERLAPS AND HAMILTONIAN ELEMENTS 
A. Overlaps 



Define 




the non- vanishing overlaps are 



{+kia \+k2(T) 



mk^a^t, |0) 



{+ka\UHB) 



z,{0\ckJl |0) 



i+ka \QE) 



^kJl |0) 
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{UHB\UHB) = {l-n%) 
{UHB\QE)^z,{l-n%) 
{QE\QE) = {l-n%) 
{-kia\-k2a) = (Olc+^Cfe^a |0) 
{-ka\LHB) = z,{0\cU^\0) 
{-ka\QH)^{0\ctS\0) 
{LHB\LHB) = n}^ 
{LHB\QH) = z,n% 
{LHB\LHB) = n% 

B. Hamiltonian Elements 



H 

Hband 
Hiocal 

Hv 



Hband + Hiocal + Hy 
ka 

r <T rgo- 



ka 
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1. H_band 



Define 



and 



XfTtr' 



j2 - 



_ v^"'ri/"i-, 

r,3.,r,5.' Jn%{l-n%)n%,{l-n%,) 
r2=riU(T'\o- ^ 



r25.,r.5.' Jn%{l - n%)n%,{l - n%,) 
ri=r2Uo-uo-' * 

ri3a,ri9<^ (1 - i^%)\/n%,{l-nj^,) 
r2=riua' J y J 



era' / J 



we will have 



V3a,V3u' /'^ 



^0 n - ) 



(+A;i(7 



band 



+ k2(T 



efcV'[5a<x'(0|Cfc,^C+^Cfc/^C+^ |0) 

k'a' 

+ (1 - (^aaO 

x(5+t(0|c,,.c+ Jt/. |0) (0|c+,,cfcv'/i/.' |0) 

+i?+,7(0|cfe,.c+ Jt/^ |0) (o|c+,,c,v'/.'/i |0) 
+S+;(0|c,,.c+ Jjt |o) (0|c+,,c,,.,/J,/., |0) 

+B;;;{o\c,,,ctJJl |o) (o|c+,,Cfc,., |o))] 
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+ (1 - 5,,,){0\CkJl |0) {xMOlK'a'^k'a'fi'fa' |0> + ^;av(0|c+^,CfeV'/a'/J. |0))] 



+(1 - S,,,){0\ckJl |0) (S++(0|c+,,Cfev'/J,/.' |0) + S+7(0|c+,,c,v'/.'& |0))] 



+(i-5,,0(i-O 

x(S+;(0|c+,,c,v'/i/.' |0) +S;,7(0|c+,,Cfev'/.'/J. |0)) 



+(l-(5,,0(l-O 

x(a;<,v(0|c+^,cfcv'/i/a' |0) + v<,v(0|c+^,Cfev//^//i |0)) 



{QE\Hband\QE) = Y.'k'ASaa'{0\fadt'aCk'afUO) 



x(B:^,{o\ct,Ak'a'flU' |0) + B:;,{o\ct,Ak'a'U'fi' m 



+(1 - 

ABt:>mtAJlh |0) (0|c+,,c,v'/i/.' |0) 

+i?:,7(0|c+,Cfe, Jt/^ |0) (0|c+,,c,v'/.'/i |0) 

+ByM^lAkJJl |0) (0|c+,,c,v'/i/.' |0) 
+B-;,{Q\clA,JJl |0) (0|c+,,c,v'/.'/J, |0»] 
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kcr Hband\LHB) = ^ek'a'{Sc7a'Za{0\CkaCt'a^k'af<j 



k'a' 



+ (1 - S„^,){0\C+J^ |0) {x^'AO\ct'a'Ck'a'fUa' |0) + V„,„{0\c+^,Ck>., |0))) 



—ka 



k'a' 



+ (1 - S,,>){0\cll |0) (i?J(0|c+,,c,v'/i/.' |0) + i?,,7(0|c+,,c,v'/.'/i |0))) 



{LHB\Hba„d\LHB) = efev'[<^aa'(0|/]a^vCfev/a |0) 

+(1 - 6,^,)n% 

x{Bt:>m,,.c,.,'flh' |0) + i^+7(0|c+,,CfeV'/.'/i |0)) 



k'a' 

+(1 - S^^>)n% 

x{Xa'a{0\c'^,^,Ck'a'fl'fa' |0) + Va'a{0\ct'^'Ck'a' fa' ft' |0)) 



{QH\Hband\QH) = J2'k'a'[Saa'{0\flct,^Ck'afa\0) 
k'a' 

+ (1 - 

x(i?,-+(0|c,V4v'/;/.' |0) + B^;,{0\ct,,,c,,,'fa'fl, |0)) 



5. HJocal 

Here we define a function for set: 

1, (7 e r 

Aa,r = { 
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a' 



Si = ^ Ernir + ^ ^ mr 

r TBa' 



r 

= y^^A^,^^yrn^m^\aS^ 



r 



r 

= y^ A^,rmr\a>S'r 



7) = y^ ■£;r(l - ^a,r) ymrmrua + y^ y^(l - ^a,r) Vmrmrug 
r 

= y^(l - ^a,r)v'mrmrutT'5'r 



r 



r 
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+ E E A'r',r[Aa,r^ - (1 - Aa,r)j^ 
= K.r;^-(l-^..r)r^lSr 



thus 



{+ha Hiocall + k2a) = {0\Ck,aC+ |0) S, + {0\Ck,Ji |0) {0\Uc+ |0) ,S25((7) 



Hlocal\UHB 



{+ka HiocailQE) = —{0\ckan\0)Ss{a) 



n 



{UHB\Hiocai\UHB)^S4{a) 



{UHB\Hiocai\QE) 
{QE\Hi,^i\QE) 



1-71°. 



n 



{-ha Hiocail - ha) = {0\cl„Ck,a |0) Si - (0|c+ J, |0) {0\flck,. |0) 525(a) 



{-ka 



Hlocal\LHB ) — 



■■mja\0)SM 



(-ka 



HlocallQH 



\ctJa\0)S,{a) 



{LHB\Hiocai\LHB) = Sj{a) 



{LHB\Hi,,,i\QH) 



n 



-S.ia) 
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{QH\Hiocai\QH) = -^Se{<j) 



3. H_V 



{+kia\Hv\ + k2(T) = ^Vk'a'[5aa'Za{{^\CkiaCl,JaCt^a |0) + (OlCfeia/lcfeVCj^^a |0)) 

fc'cr' 

+ (1 - 5.,0((0|4v/.' |0) + (0|/J,c,v' |0)) 
x(x,,,(0|c,,,c+ Jt/^ |0) + v,,,{Q\ck,.ctJJl |0))] 



{-^ka Hv\UHBj = J2Vk',,[6,.,{0\ckactJJl\0) 

k'a' 

+ (1 - 6^^,){0\ckJi |0) (w(0|c+,,A. |0) + w,,>{0\fl,Ck'a' |0))] 



k'a' 

+ + (1 - 5..0w(o|cfe./i |o) m4,,u |o) + (o|/J,CfeV' |o))] 



{UHB\Hv\UHB) = J2 yk'a'il - S.a')Vaa'{l ' O((0|c+,,A, |0) + (Ol/^Cfev |0)) 



k'a' 



{UHB\Hv\QE) = J2 Vk'a'il - Saa')(l - nJJK,,(0|c+^,/,, |0) + yaa'{0\flck'a' m 



k'a' 



{QE\Hy\QE) = ^^'-'(1 - '^-Owll - n%)m4^,U, |0) + (0|/],Cfev' |0)) 
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(-/CiCr Hv\-k2aj = ^Vk'a'[5aa'Za{{0\ct^M,JaCk2a\0) + {0\ct,jlck'aC^^^ 

k'a' 

+(1 - 5..0((o|a^v'/<x' |o) + (o|/;,Cfc,., |o)) 

{x,,,{0\cl,Ck,M, |0) + w(0|c+,Cfe, JJt |o))] 



[-ka 



Hy\LHB) = 5^VfcvKa'(0|c+Jta^v/cr|0) 
k'a' 

+(1 - 5..0(0|5t/. |0) K..(0|c+,,/., |0) + y,..(0|/;cfev' |0))] 



fc'tr' 

+ (1 - <^..')w(0|c+ |0) ((0|4,,/,, |0) + (0|/J,Cfev' |0))] 



{LHB\Hv\LHB) = - 5aa')w<((0|c+,,/,, |0) + (0|/icfcv' |0» 



k'a' 



{LHB\Hv\QH) = - <^.a')^/a(y.a'(0|4v'^' |0) + «'..'(0|/icfev' |0)) 



fe'tr' 



(g// \Hv\QH) = J] - <^aa')wn;,((0|c+,,/,, |0) + (0|/J:,Cfcv' |0)) 



k'a' 
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FIG. 1: The spectral function of electrons on the impurity site for an single orbital impurity model 
with different U and semi-circular density of states in the bath. 



FIG. 2: The spectral function of electrons on the impurity site for an single orbital impurity model 
obtained by TMA and GA lattice with U = 1. 
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FIG. 3: The density of states (DOS) obtained by DMFT+TMA for single-band Hubbard model 
on Bethe lattice at half filling. 
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FIG. 4: The density of states (DOS) obtained by DMFT+TMA of single band Hubbard model 
under U = 2 with different fillings. 
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FIG. 6: Quasi-particle weight z of single band Hubbard model obtained by DMFT+TMA versus 
U at different fillings. 
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FIG. 7: Comparison of Quasi-particle weight z for the single band Hubbard model obtained by 
DMFT+TMA, GA lattice and DMFT+ED. 
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FIG. 8: Quasi-particle weight z as the function of U for the two-band Hubbard model with SU(N) 
symmetry obtained by DMFT-I-TMA. 
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FIG. 9: Comparison of quasi-particle weight z for the two band Hubbard model with SU(N) 
symmetry at different filhngs obtained by DMFT+TMA, GA lattice and DMFT+ED. 
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FIG. 10: Comparison of quasi-particle weight z as the function of total number of particles for the 
two-band Hubbard model with SU(N) symmetry at [/ = 5 obtained by DMFT-FTMA, GA lattice 
and DMFT-hED. 
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FIG. 11: Quasi-particle weight z as the function of U for the two-degenerate-band Hubbard model 
with longitudinal Hund's coupling obtained by DMFT-I-TMA at different Jz/U . 
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FIG. 12: Comparison of quasi-particle weight z for the two-degenerate-band Hubbard model with 
longitudinal Hund's couphng obtained by DMFT-FTMA and DMFT-hED at different Jz/U. 
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CO 

FIG. 13: The spectral functions obtained by DMFT+TMA for two-nondegenerate-band Hubbard 
model with band width ratio 1 : 6 under different U with = 0.3U. 
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FIG. 14: Quasi-particle weight z of different bands as the function of U for two-nondegenerate-band 
Hubbard model with band width ratio 1 : 6 under different U with = 0.3U. 
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